Linaclotide HPLC method validation

library(readxl)

## Warning: package 'readxl' was built under R version 4.1.3

library(dplyr)

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(ggplot2)
library(plotly)

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

library(RColorBrewer)

#this.dir = rstudioapi::getActiveDocumentContext()$path
#setwd(dirname(this.dir))
linaclotide.data = read_xls(dir(pattern='XLS'), skip = 4,
                      col_types = c(
                        'text', 'skip', 'numeric', 'skip', 'text', 'text',
                        'numeric', 'skip', 'numeric', 'numeric'
                      ))

## Warning: Coercing text to numeric in H7 / R7C8: '2.442'

## Warning: Coercing text to numeric in J7 / R7C10: '172.760'

## Warning: Coercing text to numeric in K7 / R7C11: '345.369'

## Warning: Coercing text to numeric in H8 / R8C8: '2.438'

## Warning: Coercing text to numeric in J8 / R8C10: '172.931'

## Warning: Coercing text to numeric in K8 / R8C11: '345.614'

## Warning: Coercing text to numeric in H9 / R9C8: '2.440'

## Warning: Coercing text to numeric in J9 / R9C10: '172.148'

## Warning: Coercing text to numeric in K9 / R9C11: '345.825'

## Warning: Coercing text to numeric in H10 / R10C8: '2.443'

## Warning: Coercing text to numeric in J10 / R10C10: '0.602'

## Warning: Coercing text to numeric in K10 / R10C11: '1.512'

## Warning: Coercing text to numeric in H11 / R11C8: '2.442'

## Warning: Coercing text to numeric in J11 / R11C10: '110.138'

## Warning: Coercing text to numeric in K11 / R11C11: '213.756'

## Warning: Coercing text to numeric in H12 / R12C8: '2.446'

## Warning: Coercing text to numeric in J12 / R12C10: '54.947'

## Warning: Coercing text to numeric in K12 / R12C11: '102.733'

## Warning: Coercing text to numeric in H13 / R13C8: '2.450'

## Warning: Coercing text to numeric in J13 / R13C10: '12.208'

## Warning: Coercing text to numeric in K13 / R13C11: '21.988'

## Warning: Coercing text to numeric in H14 / R14C8: '2.450'

## Warning: Coercing text to numeric in J14 / R14C10: '6.031'

## Warning: Coercing text to numeric in K14 / R14C11: '10.765'

## Warning: Coercing text to numeric in H15 / R15C8: '2.451'

## Warning: Coercing text to numeric in J15 / R15C10: '1.347'

## Warning: Coercing text to numeric in K15 / R15C11: '2.403'

## Warning: Coercing text to numeric in H16 / R16C8: '2.452'

## Warning: Coercing text to numeric in J16 / R16C10: '0.660'

## Warning: Coercing text to numeric in K16 / R16C11: '1.173'

## Warning: Coercing text to numeric in H17 / R17C8: '2.453'

## Warning: Coercing text to numeric in J17 / R17C10: '0.154'

## Warning: Coercing text to numeric in K17 / R17C11: '0.267'

## Warning: Coercing text to numeric in H18 / R18C8: '2.445'

## Warning: Coercing text to numeric in J18 / R18C10: '109.991'

## Warning: Coercing text to numeric in K18 / R18C11: '213.800'

## Warning: Coercing text to numeric in H19 / R19C8: '2.447'

## Warning: Coercing text to numeric in J19 / R19C10: '57.127'

## Warning: Coercing text to numeric in K19 / R19C11: '106.836'

## Warning: Coercing text to numeric in H20 / R20C8: '2.452'

## Warning: Coercing text to numeric in J20 / R20C10: '12.077'

## Warning: Coercing text to numeric in K20 / R20C11: '21.694'

## Warning: Coercing text to numeric in H21 / R21C8: '2.452'

## Warning: Coercing text to numeric in J21 / R21C10: '6.009'

## Warning: Coercing text to numeric in K21 / R21C11: '10.728'

## Warning: Coercing text to numeric in H22 / R22C8: '2.453'

## Warning: Coercing text to numeric in J22 / R22C10: '1.245'

## Warning: Coercing text to numeric in K22 / R22C11: '2.233'

## Warning: Coercing text to numeric in H23 / R23C8: '2.455'

## Warning: Coercing text to numeric in J23 / R23C10: '0.623'

## Warning: Coercing text to numeric in K23 / R23C11: '1.108'

## Warning: Coercing text to numeric in H24 / R24C8: '2.453'

## Warning: Coercing text to numeric in J24 / R24C10: '0.135'

## Warning: Coercing text to numeric in K24 / R24C11: '0.216'

## Warning: Coercing text to numeric in H25 / R25C8: '2.446'

## Warning: Coercing text to numeric in J25 / R25C10: '101.965'

## Warning: Coercing text to numeric in K25 / R25C11: '195.745'

## Warning: Coercing text to numeric in H26 / R26C8: '2.447'

## Warning: Coercing text to numeric in J26 / R26C10: '53.896'

## Warning: Coercing text to numeric in K26 / R26C11: '100.139'

## Warning: Coercing text to numeric in H27 / R27C8: '2.451'

## Warning: Coercing text to numeric in J27 / R27C10: '11.378'

## Warning: Coercing text to numeric in K27 / R27C11: '20.314'

## Warning: Coercing text to numeric in H28 / R28C8: '2.453'

## Warning: Coercing text to numeric in J28 / R28C10: '5.725'

## Warning: Coercing text to numeric in K28 / R28C11: '10.170'

## Warning: Coercing text to numeric in H29 / R29C8: '2.453'

## Warning: Coercing text to numeric in J29 / R29C10: '1.226'

## Warning: Coercing text to numeric in K29 / R29C11: '2.200'

## Warning: Coercing text to numeric in H30 / R30C8: '2.453'

## Warning: Coercing text to numeric in J30 / R30C10: '0.589'

## Warning: Coercing text to numeric in K30 / R30C11: '1.043'

## Warning: Coercing text to numeric in H31 / R31C8: '2.453'

## Warning: Coercing text to numeric in J31 / R31C10: '0.142'

## Warning: Coercing text to numeric in K31 / R31C11: '0.246'

## Warning: Coercing text to numeric in H32 / R32C8: '2.446'

## Warning: Coercing text to numeric in J32 / R32C10: '101.758'

## Warning: Coercing text to numeric in K32 / R32C11: '194.768'

## Warning: Coercing text to numeric in H33 / R33C8: '2.446'

## Warning: Coercing text to numeric in J33 / R33C10: '107.074'

## Warning: Coercing text to numeric in K33 / R33C11: '206.443'

## Warning: Coercing text to numeric in H34 / R34C8: '2.445'

## Warning: Coercing text to numeric in J34 / R34C10: '101.591'

## Warning: Coercing text to numeric in K34 / R34C11: '194.717'

## Warning: Coercing text to numeric in H35 / R35C8: '2.452'

## Warning: Coercing text to numeric in J35 / R35C10: '11.976'

## Warning: Coercing text to numeric in K35 / R35C11: '21.640'

## Warning: Coercing text to numeric in H36 / R36C8: '2.451'

## Warning: Coercing text to numeric in J36 / R36C10: '12.022'

## Warning: Coercing text to numeric in K36 / R36C11: '21.398'

## Warning: Coercing text to numeric in H37 / R37C8: '2.451'

## Warning: Coercing text to numeric in J37 / R37C10: '11.688'

## Warning: Coercing text to numeric in K37 / R37C11: '20.788'

## Warning: Coercing text to numeric in H38 / R38C8: '2.452'

## Warning: Coercing text to numeric in J38 / R38C10: '0.462'

## Warning: Coercing text to numeric in K38 / R38C11: '0.872'

## Warning: Coercing text to numeric in H39 / R39C8: '2.453'

## Warning: Coercing text to numeric in J39 / R39C10: '0.398'

## Warning: Coercing text to numeric in K39 / R39C11: '0.711'

## Warning: Coercing text to numeric in H40 / R40C8: '2.454'

## Warning: Coercing text to numeric in J40 / R40C10: '0.388'

## Warning: Coercing text to numeric in K40 / R40C11: '0.709'

## Warning: Coercing text to numeric in H41 / R41C8: '2.455'

## Warning: Coercing text to numeric in J41 / R41C10: '0.119'

## Warning: Coercing text to numeric in K41 / R41C11: '0.220'

## Warning: Coercing text to numeric in H42 / R42C8: '2.454'

## Warning: Coercing text to numeric in J42 / R42C10: '0.120'

## Warning: Coercing text to numeric in K42 / R42C11: '0.212'

setClass('sample.prep',
         slots = list(compound='character', n.cal='numeric', reps='numeric'))
linaclotide.prep = new('sample.prep', compound='linaclotide', n.cal=7, reps=3)

n = linaclotide.prep@n.cal
lin.data.1 = linaclotide.data %>%
    slice(-c(1:4)) %>%
    cbind('stock' = c(rep(1:3, each = 7), rep(seq(4, 6), 4), 'blank')) %>%
    cbind('type' =
              c(rep('Cal', n * linaclotide.prep@reps),
                rep('QC', 12),
                'blank')) %>%
    cbind('concentration' = c(rep(c(100, 50, 10, 5, 1, 0.5, 0.1), 3),
                              rep(c(100, 10, 0.3, 0.1), each = 3),
                              0)) %>%
    rename(Area = `Peak Area`) %>%
    replace(is.na(.), 0)

lin.data.cal = lin.data.1 %>% filter(type == 'Cal')
cal.curve = lm(Area ~ concentration, data = lin.data.cal, weights = 1/concentration)
summary(cal.curve)

## 
## Call:
## lm(formula = Area ~ concentration, data = lin.data.cal, weights = 1/concentration)
## 
## Weighted Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.20650 -0.12310  0.06483  0.26783  0.59900 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    0.05765    0.06429   0.897    0.381    
## concentration  2.07752    0.01818 114.249   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.402 on 19 degrees of freedom
## Multiple R-squared:  0.9985, Adjusted R-squared:  0.9985 
## F-statistic: 1.305e+04 on 1 and 19 DF,  p-value: < 2.2e-16

cal.plot.lin = ggplot(lin.data.1, mapping = aes(concentration, Area)) +
    geom_point(aes(colour = stock, shape = type)) +
#    stat
    geom_smooth(data = lin.data.cal, 
        method = 'lm', mapping = aes(weight = 1/concentration), 
        formula = y ~ x, ) +
    labs(title='linaclotide calibration', x='concentration, ug/mL',
         y='area, mAU*s') +
    theme_classic()
#cal.plot
#
#cal.plotly = plot_ly(data = data.group, x = ~concentration, y = ~Area,
#                     color = ~stock, colors = "Reds", symbol = ~type)
#cal.plotly
cal.plot.lin

ggplotly(cal.plot.lin)

data.post.lin = lin.data.1 %>%
    group_by(stock) %>%
    mutate(calc.conc = (Area - cal.curve$coefficients[1])/
               cal.curve$coefficients[2],
            accuracy = calc.conc / concentration)
data.post.lin

data.post.lin.stock = data.post.lin %>% summarise(avg.acc = mean(accuracy))
data.post.lin.stock